Fix lost pivots in probabilistic linear algebra - #370

With the probabilistic sparse/dense linear algebra (-l 42), msolve fails on input_files/nonradical-radicalshape-31.ms for many random seeds, e.g.

msolve --random-seed 12 -t 1 -P 1 -l 42 -f input_files/nonradical-radicalshape-31.ms -o out

The output file out is empty, yet msolve exits with status 0. On macOS, 16 of seeds 1..40 fail at -t 1. -l 2 and -l 44 always succeed.

The cause is in the linear algebra. It loses pivots, so F4 computes a wrong Gröbner basis. With -v 2, the F4 run after adding the linear form shows a single round 6 x 6 ... 1 new 2 zero. This basis gives a larger degree (4 instead of 3) than the run before the linear form was added. msolve detects this inconsistency and writes no result.

Root cause

The random multipliers of the probabilistic linear algebra were drawn as

mull[j] = (int64_t)rand() & fc;

This masks with the field characteristic itself instead of reducing modulo it. For fc = 1073741827 = 2^30 + 3 (binary 1000...0011) the only possible values are 0, 1, 2, 3, 2^30, 2^30+1, 2^30+2, 2^30+3. As a result:

• the nonzero multipliers are always in {±1, ±2, ±3} mod p;
• 2^30 + 3 = fc is 0 mod p but passes the != 0 check, so about one draw in seven drops a row from the random combination.

A block stops as soon as one random combination reduces to zero, so pivots were lost, especially on small matrices. The 8- and 16-bit versions have the same problem.

Changes

• Fix random multipliers in probabilistic linear algebra: draw the multipliers as rand() % fc. Once the existing loop rejects 0, they are uniform in [1, fc-1].
• Fix buffer overflow in probabilistic_sparse_dense_echelon_form_ff_32: the multiplier array is indexed with rpb entries per thread but was allocated with ncols entries per thread. When a block has more rows than the matrix has columns, this wrote past the end of the array. It is now allocated with rpb entries per thread, as in the 8- and 16-bit versions.
• Loop over blocks in probabilistic_dense_linear_algebra_ff_{8,16,32}: the parallel loop ran up to the number of rows ntr instead of the number of blocks nb. The extra iterations had no rows to handle, so this only removes wasted work.

Tests

• New test test/diff/diff_bug-prob-sparse-dense-la.sh: runs -l 42 on nonradical-radicalshape-31 for 40 consecutive seeds with -t 1 and -t 2, and compares the results with output_files/nonradical-radicalshape-31.res. Line 5 (the random linear form, which depends on the seed) is not compared. Because it uses many seeds, the test catches the bug whatever the platform's rand(). It fails on master and passes with this PR.
• Before the fix, 16 of seeds 1..40 failed. After the fix, every run gives the same result as -l 2 for primes 1073741827, 65521 and 251, with -t 1, -t 2 and -t 4, for seeds 1..200.
• make check passes after each commit.

With this PR, the example above succeeds for all seeds: msolve writes the correct result to out and exits with status 0.

Co-Authored-By: Claude Opus 5.5 noreply@anthropic.com

Exiting with status 0 after a failed computation is a separate problem, and it is not specific to this bug. In positive characteristic, msolve does not report failures it detects, whatever their cause. Such failures can still occur after this PR: probabilistic linear algebra can still be unlucky, with a probability of about 1/p per block, and some inputs fail without genericity handling (-c 0). #371 will handle this: when msolve detects a wrong dimension or degree, it restarts the computation, and any other failure makes msolve exit with status 1.